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Abstract 

Oocytes from dairy cattle and buffaloes have severely compromised developmental competence during summer. While 
analysis of gene expression is a powerful technique for understanding the factors affecting developmental hindrance in 
oocytes, analysis by real-time reverse transcription PCR (RT-PCR) relies on the correct normalization by reference genes 
showing stable expression. Furthermore, several studies have found that genes commonly used as reference standards do 
not behave as expected depending on cell type and experimental design. Hence, it is recommended to evaluate expression 
stability of candidate reference genes for a specific experimental condition before employing them as internal controls. In 
acknowledgment of the importance of seasonal effects on oocyte gene expression, the aim of this study was to evaluate the 
stability of expression levels often well-known reference genes {ACTB, GAPDH, GUSB, HIST1H2AG, HPRTl, PPIA, RPL15, SDHA, 
TBP and YWHAZ) using oocytes collected from different categories of dairy cattle and buffaloes during winter and summer. 
A normalization factor was provided for cattle {RPL15, PPIA and GUSB) and buffaloes {YWHAZ, GUSB and GAPDH) based on 
the expression of the three most stable reference genes in each species. Normalization of non-reference target genes by 
these reference genes was shown to be considerably different from normalization by less stable reference genes, further 
highlighting the need for careful selection of internal controls. Therefore, due to the high variability of reference genes 
among experimental groups, we conclude that data normalized by internal controls can be misleading and should be 
compared to not normalized data or to data normalized by an external control in order to better interpret the biological 
relevance of gene expression analysis. 



Citation: Macabelli CH, Ferreira RM, Gimenes LU, de Carvalho NAT, Scares JG, et al. (2014) Reference Gene Selection for Gene Expression Analysis of Oocytes 
Collected from Dairy Cattle and Buffaloes during Winter and Summer. PLoS ONE 9(3): e93287. doi:10.1371/journal.pone.0093287 

Editor: Wei Yan, University of Nevada School of Medicine, United States of Amierca 

Received June 7, 2013; Accepted March 4, 2014; Published March 27, 2014 

Copyright: © 2014 Macabelli et al. This is an open-access article distributed under the terms of the Creative Commons Attribution License, which permits 
unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited. 

Funding: The research was granted by Sao Paulo Research Foundation (FAPESP; grant numbers 2009/00938-3, 2010/13384-3, 2010/09561-7, 201 1/14207-0 and 
2012/07510-1}, National Counsel of Technological and Scientific Development (CNPq; grant number 476229/201 1-1) and Coordination for the Improvement of 
Higher Level Personnel (CAPES). The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript. 

Competing interests: We would like to declare the affiliation of Marcio L. Ferraz to Vida Reprodutiva e Consultoria (Cravinhos, SP, Brazil}, of Henderson Ayres to 
MSD Saude Animal (Sao Paulo, SP, Brazil) as well as of Yeda F. Watanabe and Osnir Watanabe to WTA Tecnologia Aplicada LTDA (Cravinhos, SP, Brazil}. We confirm 
that this does not alter our adherence to all PLOS ONE policies on sharing data and materials. 

* E-mail: chiarattimr@gmail.com 

(9 These authors contributed equally to this work. 



Introduction 

Reproductive performance of bovids has long been recognized 
as a major contributor to the overall profitability of dairy 
operations. Heat stress strongly impairs female fertility in Holstein 
{Bos taurus) and buffalo {Bubalus bubalis), decreasing the herd 
reproductive efFiciency [1—7], and causing a severe economic loss, 
nearly US$ 1 billion annually [8] . It is a problem that affects about 



60% of the world cattle population [1]. Conception rates decline 
from about 40-60% in cooler months to 10-20% or lower in 
summer, depending on the severity of the thermal stress [2,9-12], 
Several studies have shown that oocytes and embryos at early 
stages of development are extremely sensitive to heat stress in both 
species [4,10,11,13-21]. Buffaloes are seasonally polyestrous with 
a marked reduction in their reproductive performance during 
summer in regions of high latitude [7,10-12,22,23]. Therefore, 
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fertility of buflFalo females is not only affected by heat stress but also 
by the sensitivity of these animals to photoperiod. However, the 
environmental effects on buffalo fertility does not appear to be due 

to a failure on follicular development, but due to an impairment of 
oocyte developmental competence [10,1 1]. The causes behind this 
cUmate-dependent developmental impairment remain widely 
unknown and studies of the effects of seasonal temperature 
oscillations on oocyte gene expression are of value in the 
understanding of the low developmental potential of bovine and 
buffalo oocytes during summer [14-16,21,24,25]. 

Gene expression analyses rely on the use of internal controls 
represented by reference genes (also known as housekeeping genes) 
that are believed to have stable expression patterns, irrespective of 
experimental treatments [26,27]. A reference gene is commonly 
used in gene expression studies to account for inherent variations 
(e.g. technical and biological) in real-time reverse transcription 
PGR (RT-PGR) studies [26,27]. Nonetheless, many experiments 
are conducted without adequate confirmation that the genes 
chosen as reference really behave as expected. Indeed, recent 
reports have indicated that genes usually chosen as internal 
controls (e.g. ACTB and GAPDH) in experiments involving somatic: 
cells do not behave as references in experiments using mammalian 
oocytes and embryos [16,28-31]. This highlights the need for 
evaluating specific candidates for each experiment, which often 
requires the selection of three or more reference genes [27-32]. 
For this purpose, the geNorm application can be applied in order 
to evaluate expression stability and determine genes, among those 
tested, that should be used [27,32]. 

Taking into consideration the importance of seasonal effects on 
oocyte gene expression, the aim of this work was to evaluate 
expression stability of ten well-known reference genes [ACTB, 
GAPDH, GUSB, HIST1H2AG, HPRTl, PPM, RPL15, SDHA, TBP 
and TWHA^ using oocytes from different categories of dairy cattle 
and buffaloes collected during winter and summer. Using the 
geNorm a normalization factor was provided for cattie (GUSB, 
PPIA and RPLIS) and buffaloes [YWHAZ, GUSB and GAPDH] 
based on expression of the three most stable reference genes in 
each species. Additionally, we provide evidence of how selection of 
internal controls is critical for gene expression analysis by 
evaluating expression of two non-reference target genes normal- 
ized by to reference genes with different expression stability. 

Materials and Methods 

All chemicals and reagents used were purchased from Sigma- 
Aldrich Ghemical Go. (St. Louis, MO, USA) unless otherwise 
stated. 

Ethics Statement 

The present study was approved by the Bioethics Gommission 
of the Faculdade de Medicina Veterinaria e Zootecnia of the 
Universidade de Sao Paulo (protocol number 1571/2008), which 
complies with the ethical principles in animal research. Gattle and 
buffaloes used in this work were provided by local farms as 
described below with consent of their owners. No specific 
permissions were required for their use because it did not involve 
endangered or protected species. We acknowledge the bovine 
farms Santa Rita - Agrindus SA (Descalvado, SP, Brazil) and Sao 
Jorge (Sao Pedro, SP, Brazil), and the Unidade de Pesquisa e 
Desenvolvimento de Registro (buffalo farm). Polo Regional do 
Vale do Ribeira-APTA (Registro, SP, Brazil) for supplying animals 
and management necessary to conduct this work. 



Sample collection 

Sampling of bovine oocytes. All bovine samples used in the 
present study were collected as reported in a previous work [33]. 
In brief, Holstein (B. taurm) cattle of three different animal 
categories were analyzed: heifers (H), peak-lactation cows (PL) and 
repeat-breeder (KB) cows. Heifers were on average 16.8±0.3 
months old, cycling, and had never been inseminated; RB were 
normal cycling lactating cows that had been inseminated several 
times (range from 4 to 1 3 services) without becoming pregnant and 
with no anatomic or infectious abnormality; and PL were normal 
cycling open cows averaging 1 10.4±3.8 days in milk. Follicle wave 
emergence was synchronized using a standard protocol [33] before 
subjecting cattle to ovum pick up (OPU). The three animal 
categories were subjected to four OPU sessions (Figure 1), two 
during winter (H, n=34; PL, n = 32; and RB, n = 31) and two 
during summer (H, n = 36; PL, n = 37; and RB, n = 36). However, 
each animal was subjected only once to OPU. A retrospective 
analysis was performed to choose the coolest and the warmest 
periods of the year in the farms where the experiments were 
conducted [33]. Gumulus-oocyte complexes (COGs) recovered 
during OPU sessions were morphologically classified based upon 
ooplasm characteristics and the number of cumulus cell layers as 
viable or unviable [33]. Approximately lO'Vn of the viable GOGs 
were randomly used for molecular analyses whereas the remaining 
90% of the COGs were used for in vitro embryo production, with 
the results published previously [33]. Oocytes used for molecular 
analyses were mechanically separated from cumulus cells by 
vortexing (3 min at maximum speed) followed by three thorough 
washes in PBS with 0.1% polyvinyl-pyrrolidone (PVP) to 
completely remove cumulus cells [16,25,30,34—38]. Expression 
of cumulus-specific genes was evaluated in denuded oocytes to 
validate the protocol used for removal of cumulus cells (FUe SI). 
Moreover, before sampling oocytes were thoroughly checked for 
the presence of somatic cells under a stereomicroscope (Figure SI). 
Oocytes were individually stored at — 80°C in 0.2 ml polystyrene 
PGR tubes with 1 ^1 of PBS with 0.1% of PVP and 1 U/^l of 
RNase inhibitor (RNase OUT, Invitrogen, Garlsbad, GA, USA). 

Sampling of buffalo oocytes. Buffalo samples used in the 
present study were collected as described above for cattle, with few 
modifications. The experiments for collection of these samples 
were conducted in a farm located in the city of Registro 
(24°30'07"S 47°50'54"W), Sao Paulo state in southeast Brazil, 
during the warm (Januar)-; summer) versus cool (July; winter) 
periods of 2011 and 2012, respectively. A retrospective analysis 
was also conducted to choose the warmest and the coolest periods 
of the year in the farm where the experiment was carried out (data 
not shown). Murrah buffaloes [B. buhalis) used in the present work 
were maintained as described in a previous report [39]. As 
performed for cattle, two categories of buffaloes were used for 
oocyte collection: nuUiparous (N) and multiparous (M). NuUipa- 
rous were normal cycling females with no anatomic or infectious 
abnormality that averaged 38.9±2.6 months old and never calved; 
multiparous cows were also cycling, calved on average 3.6±0.6 
times and were 92. 9± 10.6 months old. Multiparous buffalo cows 
were milked once daily. In order to obtain sufficient oocytes for the 
purpose of this work, foUicle wave emergence was synchronized 
using a standard protocol before OPU, similarly as performed for 
cattie [33]. This protocol is based on the use of synthetic hormones 
that stimulate the synchronized growth of follicles to antral stage 
within five days of the beginning of the treatment [40] . Briefly, on 
random days of the estrous cycle (dO), buffaloes were treated with 
2 mg of estradiol benzoate (Gonadiol, MSD Saude Animal, Gotia, 
SP, Brazil), 530 \ig of cloprostenol sodium (Giosin, MSD Saiide 
Animal) i.m. plus a norgestomet ear implant (Grestar, MSD Saiide 
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Figure 1. Schematic representation summarizing tlie experimental model. Three different animal categories of cattle [heifers (H), peak- 
lactation cows (PL) and repeat-breeder (RB)] and two of buffaloes [nulliparous (N) and multiparous (M)] were used in this study during winter (H, 
n = 34; PL, n = 32; RB, n = 31 ; N, n = 9; and M, n = 7) and summer (H, n = 36; PL, n = 37; RB, n = 36; N, n = 8; and M, n = 1 1 ). On random days of the estrous 
cycle (do), cattle were treated with 2 mg of estradiol benzoate (EB), 150 |ig of cloprostenol sodium (CS) and 50 mg of progesterone (P) whereas 
buffaloes were treated with 2 mg of EB and 530 \ig of CS. Furthermore, both cattle and buffaloes received a norgestomet ear implant for five days 
(until d5) before performing ovum pick up (OPU). Viable cumulus-oocyte complexes (COCs) recovered during OPU sessions were carefully denuded of 
cumulus cells and the oocytes individually stored for molecular analyses. 
doi:1 0.1 371/journal.pone.0093287.g001 



Animal). The norgestomet implant was maintained for five days 
(until d.l) wlien OPU was performed [40]. The two animal 
categories were subjected to two OPU sessions (Figure 1), one 
during winter (N, n = 9; and M, n = 7) and one during summer 
(N, n = 8; and M, n = 1 1). Viable COCs recovered during OPU 
sessions were carefully and thoroughly denuded of cumulus cells 
and the oocytes stored as described above for bovine samples. 

Quantification of mRNA 

RNA extraction and reverse transcription. Total RNA 

was extracted from each individual oocyte using TRIzol reagent 
(Invitrogen) according to the manufacturer's recommendations, 
with modifications. In brief, a mix containing 100 |J,1 of TRIzol 
reagent plus 1 |tl of linear acrylamide (5 Hg/M-lj Ambion Inc., 
Austin, TX) plus 9 |jl of diethylene pyrocarbonate (DEPC)-treated 
H2O (Invitrogen) was added to each sample. Samples were then 
incubated for 5 min at room temperature before adding 20 |J,1 of 
chloroform and incubating for 3 min at room temperature. This 
mixture was centrifuged at 1 2,000 x g for 15 min at 4°C and 50 [d 
of the supernatant was transferred to a new 0.2 ml polystyrene 
PGR tube to which 50 |il of isopropanol alcohol was added 
followed by incubation for 1 0 min at room temperate and 20 min 
at — 20°C. The samples were centrifuged at 12,000 x g for 10 min 
at 4°C and the supernatant removed with care (it was possible to 
visualize the pellet of linear acrylamide containing the RNA). 
Finally, ICQ nl of 75% ethanol was added into the tube, 
centrifuged at 7,500 x g for 5 min at 4°C, the supernatant 
removed (guided again by the pellet of linear acrylamide) and the 
pellet dried. 

The extracted RNA was dirc-ctly dissolved in 1 0 |il of DNase I 
solution (Invitrogen) plus 1 U/^l RNase OUT for DNA degra- 
dation as suggested by the manufacturer. To confirm the absence 
of contaminating DNA, samples were subjected to the amplifica- 
tion protocol with MT-COl primer (see protocol in [41]) before 
reverse transcription. Over one thousand copies of mitochondrial 
DNA are present in a single oocyte rendering MT-COl 
amplification an extremely sensitive method for ascertaining 
complete digestion of contaminating DNA. Finally, the RNA 
was immediately reverse transcribed into cDNA using the High- 
Capacity cDNA Reverse Transcription kit (Applied Biosystems, 



Foster City, CA, USA) according to manufacturer's protocol. 
cDNA was stored at — 20°C untU use. 

Design of TaqMan assays. A total often genes (Table SI), 
most commonly used as reference in reports of bovine oocytes 
[16,27-32,36,38,42], were selected for evaluation of bovine and 
buffalo oocytes. These genes were selected taking into account 
their roles in the cell (e.g. organization of the cytoskeleton or 
enzymatic function in energetic metabolism) in an attempt to 
select genes belonging to distinct functional classes that were not 
co-regulated. Additionally, two non-reference target genes (Table 
SI), HSPAIAB (which corresponds to amplification of both 
isoforms HSPAIA and HSPAIB) and HSP90AA1 were used to 
further understand the effect of reference gene selection on final 
gene expression data. As members of the heat shock protein 
family, expression of these two target genes is known to be affected 
by heat stress [25,43]. 

Quantification of gene transcripts was carried out using real- 
time RT-PCR using TaqMan assays. Primers and TaqMan 
probes were designed using the Primer Express software v.3.1 
(Applied Biosystems) based upon sequences available in GenBank 
(Table SI). Due to unavailability of the selected gene sequences for 
buffaloes, primers and probes used to amplify transcripts from 
cattle and buffaloes were based on bovine sequences. PGR 
products representative of each amplification assay were run onto 
2% agarose gels to check the specificity of the amplified fragment. 
Moreover, to confirm the specificity of primers to buffalo targets, 
the PGR products were sequenced on an ABI 3730 DNA Analyzer 
(using the BigDye Terminator v.3.1 Cycle Sequencing kit; Applied 
Biosystems) and the sequences aligned with the bovine sequences. 
Whenever possible, primers or probes were designed to anneal on 
exon-exon junctions and thus avoid genomic DNA amplification. 
Primers and probes were ordered from Applied Biosystems or 
Sigma- Aldrich as indicated in Table SI. 

Preamplification of target cDNAs. Before performing real- 
time RT-PCR, cDNA was preamplified using the TaqMan 
PreAmp Master Mix kit (Applied Biosystems) according to 
manufacturer's recommendations as follows: a 10-|xl reaction 
was prepared containing a pool of all primers used (each one at 
final concentration of 45 nM) plus 1 x TaqMan PreAmp Master 
Mix plus 4 nl of template (cDNA from samples), subjected to 20 
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thermal cycles and stored at — 20°C. Based on previous tests it was 
necessary to use 20 cycles of preamplification because of the low 
starting amount of cDNA regarding the target transcripts (data not 
published). The linearity of preamplification of all transcripts was 
determined as suggested by the manufacturer. Briefly, cDNA from 
pooled samples were subjected to one or two rounds of 
preamplification (20 thermal cycles each). In the latter case, 
cDNA from the first round of preamplification was diluted to 
reach the initial concentration and used as template for the second 
round of preamplification. Standard curves were generated by 
real-time RT-PCR (see below) using as template four 5-fold serial- 
dilutions of preamplified cDNAs from one or two rounds of 
preamplification. As suggested by the manufacturer, the log of 
input amount of template versus the delta (A) Ct values (Ct is the 
cycle threshold) for a serial-diluted cDNA was plotted to perform a 
finear regression. Slope values smaller than 0. 1 indicate high 
amplification efficiency, which supports the use of the 2 
method. Then, the AACt was calculated to check the linearity of 
preamplification by comparing the relative amount of a specific 
transcript after one and two rounds of amplification. As suggested 
by the manufacturer, we assumed as uniform preamplification 
AACt values smaller than ±1.5. 

Real-time reverse transcription PGR. Real-time RT-PCR 
was conducted in 15-nl reactions containing Ix TaqMan assay 
(consisting of 900 nM of primers and 250 nM of probe), 1 x 
TaqMan Gene Expression Master Mix (Applied Biosystems) and 
2 |il of template. For each sample, preamplified cDNAs were 
diluted 125 fold (1 in 125 dilution) to be used as template. Dilution 
factors were determined based on serial dilutions done to test 
amplification efficiency (see above). AH gene-specific cDNAs 
preamplified for a particular sample were run in duplicate in the 
same real-time RT-PCR plate using the ABI PRISM SDS 7500 
HT Real-Time PCR System (Applied Biosystems). The following 
cycling conditions were applied for all real-time RT-PCR 
reactions: initial denaturation at 95°C for 15 min followed by 40 
cycles consisting of 95°C for 15 sec and 60°C for 1 min. Probe 
fluorescence was read at the end of each extension step (60°C). 
Since all real-time RT-PCR assays showed high amplification 
efficiency (see above), target transcript amounts in each sample 
were linearized according to Livak and Schmittgen [26] by 2"^^^*^' 
(normahzed values) or 2 (not normalized values). The 
linearized values from averaged sample duplicates were used to 
evaluate the expression level of specific genes regarding season and 
animal category. 

GeNorm analysis 

Analysis of gene expression stability was achieved using 
qBasePLUS software v. 2.3 (Biogazelle, Zwijnaarde, Belgium) 
based on the geNorm application [32]. This approach relies on the 
principle that the expression ratio of two perfect reference genes 
should be identical in all samples, regardless of the experimental 
condition or cell type [32]. Increasing variation in this ratio 
corresponds to decreasing stability. By determining the average 
pair-wise variation between a particular candidate reference gene 
and all other candidates, the software calculates the gene stabiUty 
measure (M). Larger M values are associated with greater 
variations in gene expression [32]. By stepwise exclusion of the 
least stable gene and recalculation of the M values, the most stable 
reference genes are identified. Finally, a normalization factor is 
calculated based on the geometric mean of the expression levels of 
the best-performing reference genes [32]. To determine the 
optimal number of reference genes that should be used, the 
pairwise variation (V) was calculated between two sequential 
normalization factors (NF): NF^ and NFn+i. The NF, based on the 



geometric mean of the expression level of the n best-performing 
genes, is calculated by stepwise inclusion of an extra, less stable 
reference gene. When the gene added has a significant effect on 

the NF, a large variation is verified, which means that the gene 
should be included as reference [32]. 

Statistical analysis 

Real-time RT-PCR data were tested for normality of residuals 
and homogeneity of variances followed by ANOVA using the 
GLIMMIX procedure of SAS v. 9.3 (SAS/STAT, SAS Institute 
Inc., Gary, NC) for log normal distribution. The explanatory 
variables considered for inclusion in the models were animal 
category, season and interaction of animal category and season. 
Final results arc presented in natural log (Ln) scale (because of the 
log normal distribution considered) as normalized values of a 
specific gene transcript by the mean level of the transcript from II 
(in the case of cattie) or N (in the case of buffaloes) during winter. 
Downregulation of expression in a specific experimental group 
may be represented by negative values relative to H or N during 
winter because of Ln scale. Therefore, to avoid negative values the 
mean used for data normalization was divided by being all data 
compared against 5, which is the relative mean expression level of 
H and N during winter. Data are expressed in relation to H or N 
during winter because these are the category and season with less 
physiological (lactation) and environmental (heat stress) influence 
on fertility. Significance was considered at P<0.05. Values are 
presented as means ± standard error of the mean (S.E.M.). 

Results 

TaqMan assay screening and amplification efficiency 

Amplification specificity of all TaqMan assays used was 
confirmed by running the PCR products on agarose gels or by 
sequencing the PCR products (in the case of buffalo samples). 
Only a single assay designed to measure expression of 
HIST1H2AG failed to ampUfy its target transcript in buffalo 
samples. With this exception, targets were successfully amplified 
for all primer pairs listed in Table S 1 for both catde and buffaloes. 
Moreover, all assays used for gene expression analysis were found 
to have high efficiency of amplification. Analysis of gene 
expression used the 2 ^'^'^^ method due to the amplification 
efficiency of target and reference TaqMan assays. Efficiency was 
approximately equal for all assays tested (slope values of plotting 
amount of template versus ACt values were less than 0.1). 

Evaluation of preamplification uniformity 

Single oocytes were used for analysis of gene expression because 
of the reduced number and variable amount of oocytes available 
from donors. Samples were subjected to preamplification to 
increase the starting amount of transcripts of interest from 
individual oocytes. To check whether preamplification biased 
the relative amount of target transcripts we compared the relative 
level of all transcripts after one and two rounds of amplification by 
calculating the AACt. As a result all transcripts evaluated, except 
HSPAIAB in cattle, were found to be uniformly amplified. 
Comparison of the relative level of specific transcripts (by AACt) 
after one and two rounds of amplification were on average equal 
to 0.61 ±0.14 for bovine samples and 0.24±0.04 for buffalo 
samples (Table S2). Except for HSPAIAB in cattle that was, on 
average, equal to 2.90 ±0.43, AACt values were much lower than 
the limit 1.5, the criterion suggested by the manufacturer as 
acceptable for using preamplified cDNA in quantitative PCR. 
This corroborates the high amplification efficiency of primers used 
(see above) and justifies the use of preamplification to increase the 
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input amount of cDNA in real-time RT-PCR. In summary, these 
findings confirm that it is possible to evaluate expression of several 
genes using single oocytes by preamplifying the cDNAs of interest 
without biasing the relative amount of a transcript in relation to 
reference genes. 

Gene expression stability analysis 

Gene expression stability in oocytes from cattle and buffaloes 
was analyzed by the geNorm application (Figure 2). Of ten 
reference genes evaluated in catde, RPL15, PPIA and GUSB were 
found to be the most stable, while, in general, the remaining genes 
tested were highly variable (Figure 2A). This finding was 
confirmed by the pairwise variation analysis, as no optimal 
number of reference genes could be determined in catde 
(Figure 2a). The variability between sequential NF (based on the 
n and n+1 least variable reference targets) was relatively liigh 
(geNorm V values lower than 0.15 indicate an optimal NF). This 
means that based on the geNorm M rank (Figure 2A), even 
considering nine of ten (V''9/10) most stable reference genes in 
catde, the geNorm V index remained higher than 0.25 (Figure 2a). 
According to the geNorm instructions, because no optimal 
number of reference genes could be determined, it is recom- 
mended to use the n-1 reference genes with lowest geNorm M 
values. This means that for catde RPL15, PPIA, GUSB, ACTS, 
GAPDH, HIST1H2AG, SDIIA, TBP and HPRTl are the best 
candidates. Widi respect to buffaloes, YWHAZ GUSB and GAPDH 
were most stably expressed whereas TBP was the least stable 
(Figure 2B). Based on the pairwise variation analysis (Figure 2b), 
the optimal number of reference genes recommended in this 
experimental condition is sf'ven. As such, the optimal NF can be 
calculated as the geometric mean of YWHAZ, GUSB, GAPDH, 
ACTB, HPRTl, PPIA and RPU5. In summary, nine and seven 
reference genes, respectively, should be used for calculation of a 
suitable NF for catde and buffaloes. 

Expression profile of candidate reference genes in bovine 
and buffalo oocytes 

In order to better understand the effect of season and animal 
category' on regulation of reference genes, we investigated their 
expression levels in the oocytes from cattie and buffaloes. In cattie 
(Figure 3), most genes (^CTB, GAPDH, GUSB, HIST1H2AG, PPIA 
and RPL15] showed a significant (P<0.05) interaction of animal 
category and season, being decreased in summer compared to 
winter in H and PL, but not in RB (Figure ,3). Among these genes, 
ACTB, GAPDH and RPL15 had their expression significantiy 
decreased during winter in RB compared to H (Figure 3). 
Additionally, expression o{ ACTB, GAPDH, PPIA and RPL15 was 
increased during summer in RB compared to both H and PL. Also 
with regard to summer, expression of GUSB was increased in RB 
compared to H whereas expression of HIST1H2AG was increased 
in RB compared to PL (Figure 3). With respect to the remaining 
genes evaluated in cattie (Figure 3), expression of HPRTl, SDHA 
and TWHAZ was decreased (P<0.05) in summer compared to 
winter, regardless of animal category, whereas expression of TBP 
did not vary as a result of animal category and season. With 
respect to buffaloes (Figure 4), expression of all genes was 
significantiy decreased during summer compared to winter 
(P<0.05), but no effect of animal category or interaction was 
found. In summary, although these results give evidence of an 
overall downregulation of expression during summer, differences 
in the pattern of expression were found among the reference genes 
evaluated in cattie. 



Effect of reference gene selection on expression values of 
non-reference target genes in bovine and buffalo 
oocytes 

Aiming to understand the relevance of reference gene selection 
for gene expression experiments, we evaluated expression of two 
non-reference target genes [HSPAIAB and HSP90AA1) and 
normalized them to different NFs. With regard to cattle, 
information from HSPAIAB was excluded from analysis because 
we failed to validate it with our preamplification protocol. A 
second gene, HSP90AA1 (Figure 5), displayed an effect (P<0.05) of 
animal category and season, but not of interaction, when 
normalized by RPI15 (the most stable reference gene) or by the 
geometric mean of the three most stable references genes {RPL15, 
PPIA and GUSB). It can be inferred from tiu-sc data that 
expression of HSP90AA1 in cattie drops from winter to summer 
regardless of animal category (P<0.05). Moreover, regardless of 
season, expression of HSP90AA1 is significantiy (P<0.05) greater 
in RB than H or PL, whereas H and PL did not differ from each 
other. These results are in disagreement with those found for 
HSP90AA1 normalized by TlVHAJ^(the least stable reference gene) 
as no significant effect was detected in this case (Figure 5). 
Normalization by TBP, the only reference gene that did not vary 
among experimental groups (Figure 3), or the absence of 
normalization showed a significant (P<0.05) effect of interaction 
between animal category and season (Figure 5). This latter result 
indicates that expression of HSP90AA1 is decreased (P<0.05) in 
summer when compared to winter in H and PL, but not in RB. 
Moreover, whereas normalization by TBP resulted in greater 
(P<0.05) levels of HSP90AA1 during summer in RB when 
compared to PL, the absence of normalization resulted in greater 
(P<0.05) levels of HSP90AA1 during summer in RB when 
compared to both PL and H (Figure 5). In buffaloes, normalization 
of HSPAIAB and HSP90AA1 by the most stable reference gene 
[TWHA^ or by the geometric mean of the three most stable 
reference genes (YWHAZ, GUSB and GAPDH) yielded simflar 
results (Figure 6). Expression of HSPAIAB was significantly 
(P<0.05) greater in M than N, regardless oi s(;ason, whereas 
expression of IISP90/L41 was not affected by animal category, 
season or the interaction of them. On the other hand, 
normalization by the least stable reference gene, TBP, resulted 
in greater (P<0.05) expression of HSPAIAB in summer than 
winter, whereas no effect of animal category, season or interaction 
was found on HSP90AA1 when normalized by T5P (Figure 6). The 
absence of normalization for both HSPAIAB and HSP90AA1 
resulted in an effect (P<0.05) of season, with both decreased in 
summer compared to winter (Figure 6). In summary, these results 
provide clear evidence that reference gene selection has a major 
effect on expression values of target genes. 

Discussion 

The use of real-time RT-PCR for analysis of oocyte gene 
expression is a powerful tool for evaluating the cytoplasmic 
transcripts that may influence developmental competence. When 
considering the effect of season on oocyte viability this tool has 
been extensively employed in cattie and to a lesser extent in 
buffaloes [4,13,15,16,21,24,25]. Many of these studies in dairy 
cattie and buffaloes identified altered patterns of gene expression 
in oocytes exposed to heat stress that might be associated with 
their lower developmental competence. Studies of gene expression 
are also valuable for understanding the effect of photoperiod on 
buffalo oocytes [10,11,21]. However, analysis of gene expression 
depends on the correct selection of internal controls, as the use of 
highly variable reference genes may completely alter final gene 
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Figure 2. Gene expression stability of the candidate reference genes analyzed by the geNorm application. Average expression stability 
values (geNorm M) of evaluated genes In bovine (A) and buffalo (B) oocytes, plotted from least stable (left) to most stable (right). Pairwise variation 
analysis (geNorm V), regarding bovine (a) and buffalo (b) oocytes, between the normalization factor NFn and NF,„.,, to determine the optimal number 
of reference genes for normalization. In (a) and (b), the optimal number of reference genes is indicated by values of geNorm V below 0.1 5 (indicated 
by the dotted line). In (A), HIST1 means HIST1H2AG. 
doi:l 0.1 371/journal.pone.0093287.g002 



expression data. Thus, herein we evaluated expression stability of 
ten well-known reference genes using oocytes isolated from 
different categories of dairy cattle and buffaloes during winter 
and summer. 

Analysis of gene expression by real-time RT-PCR is the method 
of choice for studying cytoplasmic factors affecting oocyte 
competence [4,13,15,16,21,24,25,43]. Real-time RT-PCR is a 
very sensitive technique able to detect as little as a single copy of a 
transcript [44]. However, most methods available for performing 
real-time RT-PCR were developed for analysis of gene expression 
using samples recovered from a large number of .somatic cells, in 
which availability of biological material from which to extract 
RNA is not a major problem. In contrast, those working with 
oocytes and embryos usually have to deal with a limited amount of 
samples with minute amounts of RNA. As a result, real-time RT- 
PCR studies using oocytes usually do not cover many genes at a 
time and the quality of the results is low. Moreover, these studies 
usually cannot make use of several biological replicates due to the 



need for pooling samples (10-20 oocytes) to obtain rehable results 
[4,13,15,16,21,24,28-31,43]. An alternative to these problems is 
to preamplify target cDNA before carrying out real-time RT- 
PCR. Only a few studies have employed preamplification to 
analyze gene expression in oocytes or embryos [37,41,45-49], 
mosdy because of the general belief that this method biases the 
final results. However, if properly validated, preamplification is a 
valuable tool allowing reliable evaluation of several genes from 
single cells (e.g. an oocyte) [50-52]. Herein we used single oocytes 
because of the low number of oocytes recovered from individual 
donors. Moreover, because data geiicratc-d by rc-al-time RT-PCR 
are usually highly variable, the use of single oocytes enabled us to 
analyze between 10 and 20 samples per experimental group and 
thereby minimize the effect of this variation on final results. We 
show that preamplification can be accomplished using single 
oocytes from cattle and buffaloes without interfering with the 
relative amount of most transcripts. All transcripts, except 
HSPAIAB in catde, were uniformly preamphfied as the ratio 
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Figure 3. Effect of season on expression of candidate reference genes in immature oocytes from Holstein cattle. Oocytes were 
retrieved from heifers (H; n = 17 and 15 oocytes, respectively), high-producing cows In peak lactation (PL; n = 16 and 18 oocytes, respectively), and 
repeat-breeder cows (RB; n = 17 and 17 oocytes, respectively) during winter (W) and summer (S). The amounts oi ACTB (A), GAPDH (B), GUSB (C), 
HIST1H2AG (D), HPRTl (E), PPIA (F), RPL15 (G), SDHA (H), TBP (I) and YWHAZ (J) transcripts are expressed in relation to H during winter because this was 
considered as a reference group (see Material and methods for more Information). P values for animal category, season and animal category*season 
are denoted in the Insets above each graphic. Different letters over bars denote a significant difference among categories within season (P<0.05). 
*Difference between seasons within category (P<0.05). 
doi:1 0.1 371/journal.pone.0093287.g003 
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Figure 4. Effect of season on expression of candidate reference genes in immature oocytes from buffaloes. Oocytes were retrieved 
from nulliparous (N; n = 12 and 1 1 oocytes, respectively) and multlparous (M; n = 1 1 and 13 oocytes, respectively) during winter (W) and summer (S). 
The amounts oiACTB (A), GAPDH (B), GUSB (C), HPRTl (D), PPIA (E), RPL 15 (F), SDHA (G), TBP (H) and YWHAZ (I) transcripts are expressed In relation to N 
during winter because this was considered as a reference group (see Material and methods for more information). P values for animal category, 
season and animal category*season are denoted In the insets above each graphic. 
doi:1 0.1 371/journal.pone.0093287.g004 



between a specific gene and all remaining genes did not change 
after one and two rounds of preamplification. Thus, real-time RT- 
PCR could be performed with a sufficient amount of preamplified 
cDNA, enabling quantification of poorly expressed genes. 
According to the manufacturer, preamplification enables analysis 
of up to 100 genes using 1—250 ng of cDNA from single cells. This 
expands the possibilities of gene expres.sion analysis using single 
oocytes or embryos allowing for a better understanding of these 
cells without the need for pooling samples. 

Gene expression analyses are dependent on the use of internal 
controls to normalize for the number of cells used for RNA 
extraction, as well a.s taking into account degradation/loss during 
storage and extraction [27,32,53]. It is necessary to establish 
expression stability of candidate reference genes in all experimen- 
tal situations because internal controls are expected to express 
stably regardless of experimental concUtions [27,32,53]. Expression 
stability can be checked, for instance, using the geNorm 



appHcation based on the principle that the expression ratio of 
two perfect reference genes should be identical in all samples 
[27,32]. When considering internal controls for analysis of gene 
expression in oocytes and embryos several cancUdate genes should 
be tested, as expression of common reference genes can be highly 
variable in these samples [28-31,54]. Moreover, considering the 
sharp effect of season on oocyte gene expression [14— 
16,21,24,55,56], the use of multiple reference genes is recom- 
mended to measure expression levels accurately. The optimal 
number of reference genes in any specific experimental situation 
can also be calculated by the geNorm analysis [32]. Herein we 
studied ten well-known reference genes [ACTB, GAPDH, GUSB, 
HIST1H2AG, HPRTl, PPIA, RPL15, SDHA, TBP and YWHA^ 
[27-32] in oocytes retrieved from different categories of dairy 
cattle and buffsdoes during winter and summer. Based on the 
geNorm analysis, RPL15, PPIA and GUSB were the three genes 
most stable in cattle whereas YWHAZ, GUSB and GAPDH were the 
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Figure 5. Effect of reference gene selection on expression values of non-reference target genes in immature oocytes from Holstein 
cattle. Oocytes were retrieved from heifers (H; n = 1 7 and 1 5 oocytes, respectively), high-producing cows in peak lactation (PL; n = 1 6 and 1 8 oocytes, 
respectively), and repeat-breeder cows (RB; n = 17 and 17 oocytes, respectively) during winter (W) and summer (S). The amounts of HSP90AA1 (A-E) 
transcripts are expressed in relation to H during winter because this was considered as a reference group (see Material and methods for more 
information). Expression values were normalized either by YWHAZ (A), RPL15 (B), TPB (C), geometric mean of RPL15, PPIA and GUSB (D) or not 
normalized (E). P values for animal category, season and animal category*season are denoted in the insets above each graphic. Different letters over 
bars denote a significant difference among categories within season (P<0.05). *Difference between seasons within category (P<0.05). 
doi:1 0.1 371/journal.pone.0093287.g005 



three most stable in buffaloes. Curiously, YWHAZ was the most 
stable gene in buffaloes and the least stable in cattle. In general, 
the genes evaluated in this study were more stable in buffaloes than 
cattle. This is confirmed by the finding that no optimal number of 
reference genes was determined by the geNorm analysis in cattle 
whereas an optimal number of reference genes could be 
determined in buffaloes. According to this analysis, gene 
expression data from buffalo oocytes should be normalized by 
the geometric mean of YWHAZ, GUSB, GAPDH, ACTS, HPRTl, 
PPIA and RPL15. With respect to cattle, although no optimal 
number of reference genes was determined, according to the 
geNorm instructions the NF might be calculated by the geometric 
mean of nine out of ten most stable genes {RPL15, PPIA, GUSB, 
ACTB, GAPDH, HIST1H2AG, SDHA, TBP and HPRTl). The need 
for using as many as seven (buffaloes) or nine (cattle) reference 
genes for calculation of the NF gives further evidence towards the 
high variation of reference genes in these oocytes. In summary, 
these results highlight the importance of selecting appropriate 
internal controls for analysis of oocyte gene expression. 

In an attempt to compare our results with those available in the 
literature we searched for articles that similarly evaluated 
expression stability of candidate reference genes with regard to 
the effect of season or heat stress on oocyte gene expression. 
Gendelman et al. [16], evaluated expression of two genes, GAPDH 
and YWHAZ, when comparing the effect of hot and cold seasons 
on oocyte and embryo gene expression, in cattle. They found 
stable expression levels of both genes in all stages investigated, 
except at the two-cell stage in which expression of YWHAZ was 



reported to be variable. These results are in sharp contrast with 
ours as YWHAZ^ds, the least stable gene according to our analysis 
using bovine oocytes. Moreover, although GAPDH was ranked as 
the fifth most stable reference gene in our analysis, it was not very 
stable based on the geNorm M index. This result raises concern as 
GAPDH has been often used as a reference gene in experiments 
involving the effect of heat stress on oocyte and embryo gene 
expression in cattle and buffaloes [14,21,25,55]. Based on our 
results, nine and seven reference genes should be used in order to 
calculate a NF for bovine and buffalo oocytes collected during 
winter and summer. However, analysis of so many genes may be 
impractical depending on experimental conditions and goals. 
Thus, Vandesompele et al. [32] suggested the minimum use of the 
three most stable reference genes for analysis of gene expression. 
Taking this into consideration and based on our results, the NF 
might be calculated by the geometric average oi RPL15, PPIA and 
GUSB for cattle and YWHAZ, GUSB and GAPDH for buffaloes. 
Although the use of these three genes for each species does not 
meet the geNorm recommendations (geNorm V value should be 
below 0.15), their use represent a trade-off between practical 
considerations and accuracy as previously .suggested [29,32,54,57]. 
Moreover, taking into account the threshold recognized in 
software instructions (average gcNorm MS 1.0), these three genes 
chosen for each species were indeed very stable as they were all 
below the threshold. Therefore, we suggest calculating the NF 
based on geometric average of these three genes for analysis of 
oocyte gene expression from dairy cattle and buffaloes during 
winter and summer. 
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Figure 6. Effect of reference gene selection on expression values of non-reference target genes in immature oocytes from 
buffaloes. Oocytes were retrieved from nulliparous (N; n = 1 2 and 1 1 oocytes, respectively) and muciparous (M; n = 1 1 and 1 3 oocytes, respectively) 
during winter (W) and summer (S). The amounts of HSPAIA plus HSPA1B (termed in the current work as HSPAIAB; A-D) and HSP90AA1 (E-H) 
transcripts are expressed in relation to N during winter because this was considered as a reference group (see IWaterial and methods for more 
information). Expression values were normalized either by TBP (A and E), YWHAZ (B and F), geometric mean of YWHAZ, GUSB and GAPDH (C and G) or 
not normalized (D and H). P values for animal category, season and animal category*season are denoted in the insets above each graphic. 
dor:l 0.1 371 /journal.pone.0093287.g006 



Normalization by reference genes is widely employed to account 
for cell number variation in gene expression studies using somatic 
cells [27,32] where the number of cells in different samples can 
vary considerably. However, when using precise cell numbers (e.g. 
with oocytes that can be individually manipulated), there is lesser 
need for an internal control [53]. Instead, normalization can be 
done by the use of an external control (e.g. armored RNA or 
spikes) added during samphng [31,35,45,53,58]. An external 
control can account for RNA degradation and/or loss during 
sample storage and RNA extraction. Alternatively, if RNA 
degradation and loss can be minimized, even the external control 
might be omitted [41]. This represents a significant economy in 
terms of cost and time for analysis of gene expression, as three or 
more reference genes are usually employed in most experiments 
[32]. Furthermore, data that are not normalized by internal 
controls are more informative as it is possible to visualize the effect 
of the treatment on gene expression regardless of reference genes. 
This alternative is especially relevant when gene expression is 
sharply affected by treatments (e.g. effect of season). In this case, to 
avoid bias due to uneven response of internal controls to 
treatments, gene expression can be alternatively evaluated without 
normalization by reference genes. Herein we evaluated expression 
levels of reference genes without normalization by internal 
controls. Based on this analysis it was found an overall 
downregulation of expression from winter to summer in oocytes 
from catde and buffaloes. However, although the genes evaluated 
were well-known reference genes, they were differentially regulat- 
ed depending on experimental groups. Among those evaluated in 
catde, only TBP showed stable expression among experimental 
groups, HPRTl, SDHA and TWHAZhad their expression affected 
exclusively by the season and ACTS, GAPDH, GUSB, 
HIST1H2AG, PPIA and RPL15 were affected by an interaction 
of season and animal category. Even these latter genes showed 
difiFerent patterns of expression from each other, further 
highlighting an uneven regulation of reference genes in oocytes 
from cattle. In buffaloes, a more similar pattern of regulation was 
found among reference genes as the nine genes evaluated were 
expressed in similar levels within treatments. This means that for 
cattle the NF based on the geNorm analysis takes into 
consideration nine genes that expressed differentially depending 
on treatment. On the other hand, when considering the three most 
stable genes for calculation of the NF, the similar pattern of 
expression verified among them supports their use as internal 
controls. Thus, due to the uneven regulation of reference genes 
among experimental groups, data normalized by internal controls 
can be misleading and should be compared to not normalized data 
or to data normahzed by an external control to better interpret the 
biological relevance of gene expression analysis. 

Taking into account expression variability of internal controls, 
we evaluated expression of two non-reference target genes 
normalized by reference genes unevenly regulated to further 
clarify the importance of reference gene selection. For this 
purpose, HSPAIAB and HSP90AA1 were chosen as target genes 
as they are known to express differentially as a consequence of heat 
stress [25,43]. As a result, we found that different NFs can 
considerably alter expression results of HSPAIAB and HSP90AA1. 



For instance, normalization of HSP90AA1 in cattle by the least 
(TWHA^ or the most (RPL15) stable reference genes resulted in 
divergent results. If YWHA^ had been chosen as the NF, one 
would have concluded that expression of HSP90AA1 is not affected 
by animal category, season or interaction of these factors. 
Nonetheless, based on the geNorm analysis TIVHAJ^ is not 
recommended as a NF as it is highly variable. As discussed 
previously, only an effect of season was verified for TWHA^ when 
it was evaluated without normalization by internal controls. On 
the other hand, a significant interaction of animal category and 
season was found for RPL15, which may explain the differences 
found for HSP90AA1 when normalized by these two reference 
genes. Normahzation by the most stable reference gene, RPL15, or 
by the geometric mean of the three most stable reference genes 
showed that expression of HSP90AA1 is actually significantiy 
affected by animal category and season. According to these data, 
expression of HSP90AA1 drops from winter to summer regardless 
of animal category. Moreover, expression of HSP90AA1 is 
significantly greater in RB oocytes than H or PL. Yet, 
normalization by TBP or the absence of normahzation by internal 
controls indicates a significant interaction of animal category and 
season, being expression of HSP90AA1 lower in oocytes from H 
and PL during summer than winter, but not in RB. This can be 
explained by the finding that all internal controls evaluated, except 
TBP, have their expression affected by animal category and 
season. Thus, normalization by TBP or the absence of normal- 
ization highlight the overall drop in gene expression found in 
oocytes from winter to summer. On the other hand, normalization 
by TBP differed from data not normalized with regard to the 
increase oi HSP90AA1 during summer in RB when compared to H 
and PL, indicating that although statistically stable, TPB might 
have varied among experimental groups. With respect to 
buffaloes, normalization of HSP90AA1 by the least stable {TBP), 
by the most stable [YWHA^ or by the geometric mean of the three 
most stable reference genes resulted in no statistical difference. 
This is in agreement with our previous conclusion that reference 
genes beha\'ed similarly in buffaloes when analyzed without 
normalization. Moreover, compared to cattle, the geNorm analysis 
indicated that reference genes are less variable in buffaloes. In spite 
of this, a difference was noted when HSPAIAB was normalized by 
TBP or nVHAZ Normalization by TBP highlights an effect of 
season where HSPAIAB is increased in summer compared to 
winter. However, when normalized by TWHA^ this effect is no 
longer a])j)ar(;nt, while in its place, an effect of animal category 
becomes significant. This latter result indicates that HSPAIAB is 
increased in M compared to N, regardless of season, which is also 
seen when HSPAIAB is normahzed by the geometric mean of the 
three most stable reference genes. These results suggest that 
although internal controls are more stably expressed in buffaloes, 
an effect of the NF can be observed on expression of HSPAIAB. 
With respect to expression of HSPAIAB and HSP90AA1 not 
normalized by internal controls, as seen for reference genes, both 
target genes were decreased in summer compared to winter 
regardless of animal category. In summar)-, these results give clear 
evidence of an effect of the NF on expression of target genes. 
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In conclusion, a normalization factor was provided for dairy 
catde {RPL15, PPM and GUSB) and buflfaloes {YWHAZ, GUSB and 
GAPDH) based on expression of the three most stable reference 
genes in oocytes collected from each species during winter and 
summer. Normalization of non-reference target genes by these 
reference genes was shown to be considerably different compared 
to normahzation by less stable reference genes, further highUghting 
the need for careful selection of internal controls for gene 
expression analysis. 

Supporting Information 

File SI Validation of the protocol used for mechanical 
removal of cumulus cells from oocytes. Denuded oocytes 
used in molecular analyses were evaluated regarding the presence 
of cumulus-specific transcripts to confirm the absence of 
contaminating cumulus cells. 
(DOC) 

Figure SI Comparison of chemical and mechanical 
removal of cumulus cells. The pictures depict oocytes that 
had the zona pellucida (ZP) chemically removed by treatment with 
0.1% (w/v) of pronase for 5 min (A) and oocytes that were 
mechanically separated from cumulus cell by vortexing for 3 min 
at maximum speed (B). Note that no cumulus cell is attached to 
oocyte surface regardless of the presence of the ZP. Images were 
taken at 400 x magnification. 
(TIF) 

Figure S2 Effect of chemical and mechanical removal of 
cumulus cells on expression of cumulus-specific genes. 

Cumulus-oocyte complexes (COCs) were denuded of cumulus 
cells by chemical treatment to remove the zona pellucida (without 
ZP) or by vortexing (with ZP). The amounts of KITLG, EGFR and 
FSHR transcripts are expressed in relation to ZP-free oocytes. 
Expression values were normalized by the geometric mean of 
RPL15, PPIA and GUSB. P values refer to comparisons made for 
each gene between oocytes with and without the ZP. Transcripts 
from FSHR were not found in both groups. 
(TIF) 

Figure S3 Effect of chemical and mechanical removal of 
cumulus cells on expression of cumulus-specific and 
reference genes. Cumulus-oocyte complexes (COCs) were 
denuded of cumulus cells by chemical treatment to remove the 
zona pellucida (without ZP) or by vortexing (with ZP). The 
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